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o. ABSTRACT 

Q^. We present six simulations of Galactic stellar haloes formed by the tidal disruption of 

accreted dwarf galaxies in a fully cosmological setting. Our model is based on the Aquarius 
project, a suite of high resolution N-body simulations of individual dark matter haloes. We 
^ tag subsets of particles in these simulations with stellar populations predicted by the GAL- 

^ I FORM semi-analytic model. Our method self-consistently tracks the dynamical evolution and 

disruption of satellites from high redshift. The luminosity function and structural properties of 
surviving satellites, which agree well with observations, suggest that this technique is appro- 
priate. We find that accreted stellar haloes are assembled between 1 < z < 7 from less than 
5 significant progenitors. These progenitors are old, metal-rich satellites with stellar masses 
similar to the brightest Milky Way dwarf spheroidals (10^ — 10^ Mq). In contrast to previous 
\ stellar halo simulations, we find that several of these major contributors survive as self-bound 

. systems to the present day. Both the number of these significant progenitors and their in- 

fall times are inherently stochastic. This results in great diversity among our stellar haloes, 
which amplifies small differences between the formation histories of their dark halo hosts. 
The masses 10^ — 10^ Mq) and density/surface-brightness profiles of the stellar haloes 
(from 10-100 kpc) are consistent with expectations from the Milky Way and M31. Each 
halo has a complex structure, consisting of well-mixed components, tidal streams, shells and 
other subcomponents. This structure is not adequately described by smooth models. The cen- 
KJ{ \ tral regions (< 10 kpc) of our haloes are highly prolate (c/a ~ 0.3), although we find one 

. example of a massive accreted thick disc. Metallicity gradients in our haloes are typically sig- 

■ nificant only where the halo is built from a small number of satellites. We contrast the ages 

and metallicities of halo stars with surviving satellites, finding broad agreement with recent 
observations. 
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1 INTRODUCTION 

An extended and diffuse stellar halo envelops the Milky Way. Al- 
though only an extremely small fraction of the stars in the Solar 
neighbourhood belong to this halo, they can be easily recognized 
by their extreme kinematics and metallicities. Stellar populations 
with these properties can now be followed to distances in excess of 
100 kpc using luminous tracers such as RR Lyraes, blue horizon- 
tal branch stars, metal-poor giants and globular clusters (e.g. Oort 
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1926; Baade 1944; Eggen, Lynden-Bell & Sandage 1962; Searle & 
Zinn 1978; Yanny et al. 2000; Vivas & Zinn 2006; Morrison et al. 
2009). 

In recent years, large samples of halo-star velocities (e.g. Mor- 
rison et al. 2000; Starkenburg et al. 2009) and 'tomographic' photo- 
metric and spectroscopic surveys have shown that the stellar halo is 
not a single smoothly-distributed entity, but instead a superposition 
of many components (Belokurov et al. 2006; Juric et al. 2008; Bell 
et al. 2008; Carollo et al. 2007, 2009; Yanny et al. 2009). Notable 
substructures in the Milky Way halo include the broad stream of 
stars from the disrupting Sagittarius dwarf galaxy (Ibata, Gilmore 
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& Irwin 1994; Ibata et al. 2001), extensive and diffuse overden- 
sities (Juric et al. 2008; Belokurov et al. 2007a; Watkins et al. 
2009), the Monoceros 'ring' (Newberg et al. 2002; Ibata et al. 2003; 
Yanny et al. 2003), the orphan stream (Belokurov et al. 2007b) and 
other kinematically cold debris (Schlaufman et al. 2009). Many of 
these features remain unclear. At least two kinematically distinct 
'smooth' halo components have been identified from the motions of 
stars in the Solar neighbourhood, in addition to one or more 'thick 
disc' components (Carollo et al. 2009). Although current obser- 
vations only hint at the gross properties of the halo and its sub- 
structures, some general properties are well-established: the halo is 
extensive (> 100 kpc), metal-poor (([Fe/H]) ^ —1.6, e.g. Laird 
et al. 1988; Carollo et al. 2009) and contains of the order of 0.1-1% 
of the total stellar mass of the Milky Way (recent reviews include 
Freeman & Bland-Hawthorn 2002; Helmi 2008). 

Low surface-brightness features seen in projection around 
other galaxies aid in the interpretation of the Milky Way's stellar 
halo, and vice versa. Diffuse concentric 'shells' of stars on 100 kpc 
scales around otherwise regular elliptical galaxies have been at- 
tributed to accretion events (e.g. Schweizer 1980; Quinn 1984). 
Recent surveys of M31 (e.g. Ferguson et al. 2002; Kalirai et al. 
2006; Ibata et al. 2007; McConnachie et al. 2009) have revealed 
an extensive halo (to ^ 150 kpc) also displaying abundant sub- 
structure. The surroundings of other nearby Milky Way analogues 
are now being targeted by observations using resolved star counts 
to reach very low effective surface brightness limits, although as 
yet no systematic survey has been carried out to sufficient depth, 
(e.g. Zibetti & Ferguson 2004; McConnachie et al. 2006; de Jong, 
Radburn-Smith & Sick 2008; Barker et al. 2009; Ibata, Mouhcine, 
& Rejkuba 2009). A handful of deep observations beyond the Lo- 
cal Group suggest that stellar haloes are ubiquitous and diverse 
(e.g. Sackett et al. 1994; Shang et al. 1998; Malin & Hadley 1999; 
Martmez-Delgado et al. 2008, 2009; Faundez-Abans et al. 2009). 

Stellar haloes formed from the debris of disrupted satellites are 
a natural byproduct of hierarchical galaxy formation in the ACDM 
cosmolog30. The entire assembly history of a galaxy may be en- 
coded in the kinematics, metallicities, ages and spatial distributions 
of its halo stars. Even though these stars constitute a very small 
fraction of the total stellar mass, the prospects are good for recov- 
ering such information from the haloes of the Milky Way, M3 1 and 
even galaxies beyond the Local Group (e.g. Johnston, Hernquist 
& Bolte 1996; Helmi & White 1999). In this context, theoretical 
models can provide useful 'blueprints' for interpreting the great di- 
versity of stellar haloes and their various sub-components, and for 
relating these components to fundamental properties of galaxy for- 
mation models. Alongside idealised models of tidal disruption, ab 
initio stellar halo simulations in realistic cosmological settings are 
essential for direct comparison with observational data. 

In principle, hydrodynamical simulations are well- suited to 
this task, as they incorporate the dynamics of a baryonic compo- 
nent self-consistently. However, many uncertainties remain in how 
physical processes such as star formation and supernova feedback, 
which act below the scale of individual particles or cells, are im- 
plemented in these simulations. The computational cost of a sin- 

^ In addition to forming components of the accreted stellar halo, infalling 
satellites may cause dynamical heating of a thin disc formed 'in situ' 
(e.g. Toth & Ostriker 1992; Velazquez & White 1999; Benson et al. 2004; 
Kazantzidis et al. 2008) and may also contribute material to an accreted 
thick disc (Abadi, Navarro & Steinmetz 2006) or central bulge. We discuss 
these additional contributions to the halo, some of which are not included 
in our modelling, in Section [33] 



gle state-of-the-art hydrodynamical simulation is extremely high. 
This cost severely limits the number of simulations that can be 
performed, restricting the freedom to explore different parameter 
choices or alternative assumptions within a particular model. The 
computational demands of hydrodynamical simulations are com- 
pounded in the case of stellar halo models, in which the stars of 
interest constitute only 1% of the total stellar mass of a Milky 
Way-like galaxy. Even resolving the accreted dwarf satellites in 
which a significant proportion of these halo stars may originate is 
close to the limit of current simulations of disc galaxy formation. To 
date, few hydrodynamical simulations have focused explicitly on 
the accreted stellar halo (recent examples include Bekki & Chiba 
2001; Brook et al. 2004; Abadi et al. 2006 and Zolotov et al. 2009). 

In the wider context of simulating the 'universal' popula- 
tion of galaxies in representative (> 100 Mpc^) cosmological vol- 
umes, these practical limitations of hydrodynamical simulations 
have motivated the development of a powerful and highly success- 
ful alternative, which combines two distinct modelling techniques: 
well-understood high-resolution N-body simulations of large-scale 
structure evolved self-consistently from ACDM initial conditions 
and fast, adaptable semi-analytic models of galaxy formation with 
very low computational cost per run (Kauffmann, Nusser & Stein- 
metz 1997; Kauffmann et al. 1999; Springel et al. 2001; Hatton 
et al. 2003; Kang et al. 2005; Springel et al. 2005; Bower et al. 
2006; Croton et al. 2006; De Lucia et al. 2006). In this paper we de- 
scribe a technique motivated by this approach which exploits com- 
putationally expensive, ultra-high-resolution N-body simulations 
of individual dark matter haloes by combining them with semi- 
analytic models of galaxy formation. Since our aim is to study the 
spatial and kinematic properties of stellar haloes formed through 
the tidal disruption of satellite galaxies, our technique goes beyond 
standard semi- analytic treatments. 

The key feature of the method presented here is the dynamical 
association of stellar populations (predicted by the semi- analytic 
component of the model) with sets of individual particles in the 
N-body component. We will refer to this technique as 'particle tag- 
ging'. We show how it can be applied by combining the Aquarius 
suite of six high resolution isolated ^10^^ Mq dark matter haloes 
(Springel et al. 2008a,b) with the GALFORM semi-analytic model 
(Cole et al. 1994, 2000; Bower et al. 2006). These simulations can 
resolve structures down to ^ 10^ Mq, comparable to the least mas- 
sive dark halo hosts inferred for Milky Way satellites (e.g. Strigari 
et al. 2007; Walker et al. 2009). 

Previous implementations of the particle-tagging approach 
(White & Springel 2000; Diemand, Madau & Moore 2005; Moore 
et al. 2006; Bullock & Johnston 2005; De Lucia & Helmi 2008) 
have so far relied on cosmological simulations severely limited by 
resolution (Diemand et al. 2005; De Lucia & Helmi 2008) or else 
simplified higher resolution N-body models (Bullock & Johnston 
2005). In the present paper, we apply this technique as a postpro- 
cessing operation to a 'fully cosmological' simulation, in which 
structures have grown ab initio, interacting with one another self- 
consistently. The resolution of our simulations is sufficient to re- 
solve stellar halo substructure in considerable detail. 

With the aim of presenting our modelling approach and ex- 
ploring some of the principal features of our simulated stellar 
haloes, we proceed as follows. In Section[2]we review the Aquarius 
simulations and their post-processing with the GALFORM model, 
and in Section [3] we describe our method for recovering the spa- 
tial distribution of stellar populations in the halo by tagging parti- 
cles. We calibrate our model by comparing the statistical properties 
of the surviving satellite population to observations; the focus of 



Galactic stellar haloes in the CDM model 3 



this paper is on the stellar halo, rather on than the properties of 
these satellites. In Section |4l we describe our model stellar haloes 
and compare their structural properties to observations of the Milky 
Way and M31. We also examine the assembly history of the stellar 
haloes in detail (Section I4l2l ) and explore the relationship between 
the haloes and the surviving satellite population. Finally, we sum- 
marise our results in Section[5] 



2 AQUARIUS AND GALFORM 

Our model has two key components: the Aquarius suite of six 
high-resolution N-body simulations of Milky Way-like dark mat- 
ter haloes, and GALFORM, a semi-analytic model of galaxy forma- 
tion. The technique of post-processing an N-body simulation with 
a semi-analytic model is well established (Kauffmann et al. 1999; 
Springel et al. 2001; Helly et al. 2003; Hatton et al. 2003; Kang 
et al. 2005; Bower et al. 2006; De Lucia et al. 2006), although 
its application to high-resolution simulations of individual haloes 
such as Aquarius is novel and we review relevant aspects of the 
GALFORM code in this context below. 

Here, in the post-processing of the N-body simulation, the 
stellar populations predicted by GALFORM to form in each halo 
are also associated with 'tagged' subsets of dark matter particles. 
By following these tagged dark matter particles, we track the evolv- 
ing spatial distribution and kinematics of their associated stars, in 
particular those that are stripped from satellites to build the stel- 
lar halo. This level of detail regarding the distribution of halo stars 
is unavailable to a standard semi-analytic approach, in which the 
structure of each galaxy is represented by a combination of ana- 
lytic density profiles. 

Tagging particles in this way requires the fundamental as- 
sumption that baryonic mass nowhere dominates the potential and 
hence does not perturb the collisionless dynamics of the dark mat- 
ter. Generally, a massive thin disc is expected to form at some point 
in the history of our 'main' haloes. Although our semi-analytic 
model accounts for this thin disc consistently, our dark matter tag- 
ging scheme cannot represent its dynamics. For this reason, and 
also to avoid confusion with our accreted halo stars, we do not at- 
tempt to tag dark matter to represent stars forming in situ in a thin 
disc at the centre of the main halo. The approximation that the dy- 
namics of stars can be fairly represented by tagging dark matter par- 
ticles is justifiable for systems with high mass-to-light ratios such 
as the dwarf satellites of the Milky Way and M31 (e.g. Simon & 
Geha 2007; Walker et al. 2009), the units from which stellar haloes 
are assembled in our models. 

2.1 The Aquarius Haloes 

Aquarius (Springel et al. 2008a) is a suite of high-resolution sim- 
ulations of six dark matter haloes having masses within the range 
1 — 2 X 10^^ Mq, comparable to values typically inferred for the 
Milky Way halo (e.g. BattagUa et al. 2005; Smith et al. 2007; Li 
& White 2008; Xue et al. 2008). By matching the abundance of 
dark matter haloes in the Millennium simulation to the SDSS stel- 
lar mass function, Guo et al. (2009) find 2.0 x 10^^ Mq (with a 
10-90% range of 0.8 x IO^^Mq to 4.7 x IO^^Mq). This value is 
sensitive to the assumption that the Milky Way is a typical galaxy, 
and to the adopted Milky Way stellar mass (5.5 X 10^° M0;Flynn 
et al. 2006). 

The Aquarius haloes were selected from a lower resolution 
version of the Millennium-II simulation (Boylan-Kolchin et al. 



Table 1. Properties of the six Aquarius dark matter halo simulations 
(Springel et al. 2008a) on which the models in this paper are based. The 
first column labels the simulation (abbreviated from the as Aq-A-2, Aq-B- 
2 etc.). From left to right, the remaining columns give the particle mass 
mp, the number of particles within r2oo, the virial radius at 2; = 0; the 
virial mass of the halo, M200; and the maximum circular velocity, Vmax, 
and corresponding radius, rmax- Virial radii are defined as the radius of a 
sphere with mean inner density equal to 200 times the critical density for 
closure. 
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2009) and individually resimulated using a multi-mass particle 
('zoom') technique. In this paper we use the 'level 2' Aquar- 
ius simulations, the highest level at which all six haloes were 
simulated. We refer the reader to Springel et al. (2008a,b) for a 
comprehensive account of the entire simulation suite and demon- 
strations of numerical convergence. We list relevant properties of 
each halo/simulation in Table [T] The simulations were carried 
out with the parallel Tree-PM code GADGET- 3, an updated ver- 
sion of GADGET- 2 (Springel 2005). The Aq-2 simulations used a 
fixed comoving Plummer-equivalent gravitational softening length 
of e = 48 pc. ACDM cosmological parameters were adopted 
as Qm = 0.25, Qa — 0.75, ag — 0.9, ris — 1, and Hubble con- 
stant Hq — lOO/i kms~^Mpc~^. A value ofh = 0.73 is assumed 
throughout this paper. These parameters are identical to those used 
in the Millennium Simulation and are marginally consistent with 
WMAP 1- and 5-year constraints (Spergel et al. 2003; Komatsu 
et al. 2009). 



2.2 The GALFORM Model 

N-body simulations of cosmic structure formation supply infor- 
mation on the growth of dark matter haloes, which can serve as 
the starting point for a semi- analytic treatment of baryon accre- 
tion, cooling and star formation (see Baugh 2006, for a compre- 
hensive discussion of the fundamental principles of semi-analytic 
modelling). The Durham semi-analytic model, GALFORM, is used 
in this paper to postprocess the Aquarius N-body simulations. The 
GALFORM code is controlled by a number of interdependent pa- 
rameters which are constrained in part by theoretical limits and re- 
sults from hydrodynamical simulations. Remaining parameter val- 
ues are chosen such that the model satisfies statistical comparisons 
with several datasets, for example the galaxy luminosity function 
measured in several wavebands (e.g. Baugh et al. 2005; Bower et al. 
2006; Font et al. 2008). Such statistical constraints on large scales 
do not guarantee that the same model will provide a good descrip- 
tion of the evolution of a single 'Milky Way' halo and its satellites. 
A model producing a satellite galaxy luminosity function consis- 
tent with observations is a fundamental prerequisite for the work 
presented here, in which a proportion of the total satellite popula- 
tion provides the raw material for the assembly of stellar haloes. We 
demonstrate below that the key processes driving galaxy formation 
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on small scales are captured to good approximation by the existing 
GALFORM model and parameter values of Bower et al. (2006). 

Many of the physical processes of greatest relevance to galaxy 
formation on small scales were explored within the context of semi- 
analytic modelling by Benson et al. (2002b). Of particular signifi- 
cance are the suppression of baryon accretion and cooling in low 
mass haloes as the result of photoheating by a cosmic ionizing 
background, and the effect of supernova feedback in shallow po- 
tential wells. Together, these effects constitute a straightforward 
astrophysical explanation for the disparity between the number of 
low mass dark subhaloes found in N-body simulations of Milky 
Way-mass hosts and the far smaller number of luminous satellites 
observed around the Milky Way (the so-called 'missing satellite' 
problem). Recent discoveries of faint dwarf satellites and an im- 
proved understanding of the completeness of the Milky Way sam- 
ple (Koposov et al. 2008; Tollerud et al. 2008, and refs. therein) 
have reduced the deficit of observed satellites, to the point of qual- 
itative agreement with the prediction of the model of Benson et al. 
(2002b). At issue now is the quality (rather than the lack) of agree- 
ment between such models and the data. We pay particular attention 
to the suppressive effect of photoheating. This is a significant pro- 
cess for shaping the faint end of the satellite luminosity function 
when, as we assume here, the strength of supernova feedback is 
fixed by constraints on the galaxy population as a whole. 

2.2.1 Reionization and the satellite luminosity function 

A simple model of reionization heating based on a halo mass de- 
pendent cooling threshold (Benson et al. 2003) is implemented in 
the Bower et al. (2006) model of GALFORM. This threshold is set 
by parameters termed Vent and Zcnt- No gas is allowed to cool 
within haloes having a circular velocity below Kut at redshifts be- 
low ZcMt- To good approximation, this scheme reproduces the link 
between the suppression of cooling and the evolution of the 'filter- 
ing mass' (as defined by Gnedin 2000) found in the more detailed 
model of Benson et al. (2002b), where photoheating of the inter- 
galactic medium was modelled explicitly. In practice, in this sim- 
ple model, the value of Kut is most important. Variations in ^cut 
within plausible bounds have a less significant effect on the 2; = 
luminosity function. 

As stated above, we adopt as a fiducial model the GALFORM 
implementation and parameters of Bower et al. (2006). However, 
we make a single parameter change, lowering the value of Kut 
from 50kms~^ to 30kms~^. This choice is motivated by re- 
cent ab initio cosmological galaxy formation simulations incorpo- 
rating the effects of photoionization self-consistently (Hoeft et al. 
2006; Okamoto, Gao & Theuns 2008; Okamoto & Frenk 2009; 
Okamoto et al. 2009). These studies find that values of Kut ^ 
25 — 35kms~^ are preferable to the higher value suggested by 
the results of Gnedin (2000) and adopted in previous semi-analytic 
models (e.g. Somerville 2002; Bower et al. 2006; Croton et al. 
2006; Li, De Lucia & Helmi 2009a). Altering this value affects 
only the very faint end of the galaxy luminosity function, and so 
does not change the results of (Bower et al. 2006). The choice of a 
fiducial set of semi-analytic parameters in this paper illustrates the 
flexibility of our approach to modelling stellar haloes. The N-body 
component of our models - Aquarius - represents a considerable 
investment of computational time. In contrast, the semi-analytic 
post-processing can be re-run in only a few hours, and can be easily 
'upgraded' (by adding physical processes and constraints) in order 
to provide more detailed output, explore the consequences of pa- 
rameter variations, or compare alternative semi-analytic models. 




Figure 1. The cumulative V-band luminosity functions (LPs) of satellite 
galaxies for the six Aquarius haloes, adopting in GALFORM the parameters 
of Bower et al. (2006) with Vent = 30 kms~^. These LPs include the ef- 
fects of tidal stripping measured from our assignment of stars to dark matter 
particles (Section [3), although this makes only a small difference to the LP 
from our semi-analytic model alone. All galaxies within 280 kpc of the halo 
centre are counted as satellites (the total number of contributing satellites 
in each halo is indicated in the legend). The stepped line (grey, with error 
bars) shows the observed mean luminosity function found by Koposov et al. 
(2008) for the MW and M31 satellite system (also to 280 kpc), assuming 
an NPW distribution for satellites in correcting for SDSS sky coverage and 
detection efficiency below Mv = — 10. The colour-coding of our haloes in 
this figure is used throughout. 



The V-band satellite luminosity function resulting from the 
application of the GALFORM model described above to each Aquar- 
ius halo is shown in Fig. [T] Satellites are defined as all galax- 
ies within a radius of 280 kpc from the centre of potential in the 
principal halo, equivalent to the limiting distance of the Koposov 
et al. (2008) completeness-corrected observational sample. These 
luminosity functions are measured from the particle realisations 
of satellites that we describe in the following section, and not di- 
rectly from the semi- analytic model. They therefore account for 
the effects of tidal stripping, although these are minor: the fraction 
of satellites brighter than My = —10 is reduced very slightly in 
some of the haloes. In agreement with the findings of Benson et al. 
(2002a), the model matches the faint end of the luminosity function 
well, but fewer bright satellites are found in each of our six models 
than are observed in the mean of the Milky Way -1- M31 system, 
although the number of objects concerned is small. The true abun- 
dance of bright satellites for Milky Way-mass hosts is poorly con- 
strained at present, so it is unclear whether or not this discrepancy 
reflects cosmic variance, a disparity in mass between the Aquarius 
haloes and the Milky Way halo, or a shortcoming of our fiducial 
Bower et al. (2006) model. A modification of this model in which 
the tidal stripping of hot gas coronae around infalling satellites is 
explicitly calculated (rather than assuming instantaneous removal; 
see Font et al. 2008) produces an acceptable abundance of bright 
satellites. 
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2. 2. 2 Further details 

Within GALFORM, cold gas is transferred from tidally destroyed 
satellites to the disc of the central galaxy when their host subhaloes 
are no longer identified at the resolution limit imposed by SUB- 
FIND. In the Aq-2 simulations this corresponds to a minimum re- 
solved dark halo mass of ^ 3 x IO^Mq. In the GALFORM model 
of Bower et al. (2006), which does not include tidal stripping or 
a 'stellar halo' component, the satellite galaxy is considered to be 
fully disrupted (merged) at this point: its stars are transferred to 
the bulge component of the central galaxy. By contrast, our parti- 
cle representation (described in Section [3]) allows us to follow the 
actual fate of the satellite stars independently of this choice in the 
semi-analytic model. This choice is therefore largely a matter of 
'book-keeping'; we have ensured that adopting this approach does 
not prematurely merge galaxies in the semi-analytic model that are 
still capable of seeding new stellar populations into the particle rep- 
resentation. Semi-analytic models based on N-body simulations of- 
ten choose to 'follow' satellites with dark haloes falling below the 
numerical resolution by calculating an appropriate merger time- 
scale from the last-known N-body orbital parameters, accounting 
for dynamical friction. However, the resolution of Aquarius is suf- 
ficiently high to make a simpler and more self-consistent approach 
preferable in this case, preserving the one-to-one correspondence 
between star-forming semi-analytic galaxies and bound objects in 
the simulation. We have checked that allowing semi-analytic galax- 
ies to survive without resolved subhaloes, subject to the treatment 
of dynamical friction used by Bower et al. (2006), affects only the 
faintest (Mv ^ 0) part of the survivor luminosity function. The true 
nature and survival of these extremely faint sub-resolution galaxies 
remains an interesting issue to be addressed by future semi-analytic 
models of galactic satellites. 

In Table [2 (Section 131) we list the V-band magnitudes and total 
stellar masses of the central galaxies that form in the six Aquar- 
ius haloes. A wide range is evident, from an M31 -analogue in halo 
Aq-C, to an M3 3 -analogue in Aq-E. This is not unexpected: the 
Aquarius dark haloes were selected only on their mass and isola- 
tion, and these criteria alone do not guarantee that they will host 
close analogues of the Milky Way. The scaling and scatter in the 
predicted relationship between halo mass and central galaxy stellar 
mass are model-dependent. With the GALFORM parameter values 
of Bower et al., the mean central stellar mass in a typical Aquarius 
halo (Mhaio 1.4 X IO^^Mq) is - 1.5 x 10^° Mq, approxi- 
mately a factor of 3-4 below typical estimates of the stellar mass 
of the Milky Way (- 6 x 10^° Mq Flynn et al. 2006); the scat- 
ter in Mgai for our central galaxies reflects the overall distribution 
produced by the model of Bower et al. (2006) for haloes of this 
mass. The model of De Lucia et al. (2006), which like the Bower 
et al. (2006) model was constrained using statistical properties of 
bright field and cluster populations, produces a mean central stellar 
mass of - 4 X 10^° Mq for the typical halo mass of the Aquarius 
simulations, as well as a smaller scatter about the mean value. 

In light of these modelling uncertainties and observational un- 
certainties in the determination of the true Milky Way dark halo 
mass to this precision, we choose not to scale the Aquarius haloes 
to a specific mass for 'direct' comparison with the Milky Way. The 
results we present concerning the assembly and structure of stel- 
lar haloes and the ensemble properties of satellite systems should 
not be sensitive to whether or not their galaxies are predicted to 
be direct analogues of the Milky Way by the Bower et al. (2006) 
GALFORM model. Therefore, in interpreting the absolute values of 
quantities compared to observational data in the following sections. 



it should be borne in mind that we model a range of halo masses 
that could lie somewhat below the likely Milky Way value. 

The Bower et al. (2006) implementation of GALFORM results 
in a mass-metallicity relation for faint galaxies which is slightly 
steeper than that derived from the satellites of the Milky Way and 
M31 (e.g. Mateo 1998; Kirby et al. 2008; see also Tremonti et al. 
2004 and refs. therein). This results in model galaxies being on 
average 0.5 dex more metal-poor in [Fe/H] than the observed 
relation at magnitudes fainter than My ^ —10. Whilst it would 
be straightforward to make ad hoc adjustments to the model pa- 
rameters in order to match this relation, doing so would violate 
the agreement established between the Bower et al. (2006) param- 
eter set and a wide range of statistical constraints from the bright 
(Mv < — 19) galaxy population. 



3 BUILDING STELLAR HALOES 
3.1 Assigning Stars To Dark Matter 

Observations of the stellar velocity distributions of dwarf 
spheroidal satellites of the Milky Way imply that these objects are 
dispersion- supported systems with extremely high mass-to-light ra- 
tios, of order 10-1000 (e.g. Mateo 1998; Simon & Geha 2007; Stri- 
gari et al. 2007; Wolf et al. 2009; Walker et al. 2009). As we de- 
scribe in this section, in order to construct basic models of these 
high-M/L systems without simulating their baryon dynamics ex- 
plicitly, we will assume that their stars are formed 'dynamically 
coupled' to a strongly bound fraction of their dominant dark mat- 
ter component, and will continue to trace that component through- 
out the simulation. Here we further assume that the depth at which 
stars form in a halo potential well depends only on the total mass 
of the halo. While these assumptions are too simplistic a descrip- 
tion of stellar dynamics in such systems to compare with detailed 
structural and kinematic observations, we show that they none the 
less result in half-light radii and line-of-sight velocity dispersions in 
agreement with those of Milky Way dwarf spheroidals. Hence the 
disruption of a fraction of these model satellites by tidal forces in 
the main halo should reproduce stellar halo components ('streams') 
at a level of detail sufficient for an investigation of the assembly and 
gross structure of stellar haloes. We stress that these comparisons 
are used as constraints on the single additional free parameter in 
our model, and are not intended as predictions of a model for the 
satellite population. 

In the context of our GALFORM model, the stellar content of a 
single galaxy can be thought of as a superposition of many distinct 
stellar populations, each defined by a particular formation time and 
metallicity. Although the halo merger tree used as input to GAL- 
FORM is discretized by the finite number of simulation outputs 
(snapshots), much finer interpolating timesteps are taken between 
snapshots when solving the differential equations governing star 
formation. Consequently, a large number of distinct populations 
are 'resolved' by GALFORM. However, we can update our parti- 
cle (dynamical) data (and hence, can assign stars to dark matter) 
only at output times of the pre-existing N-body simulation. For the 
purposes of performing star-to-dark-matter assignments we reduce 
the fine-grained information computed by GALFORM between one 
output time and the next to a single aggregated population of 'new 
stars' formed at each snapshot. 

As discussed above and in Section [T] we adopt the fundamen- 
tal assumption that the motions of stars can be represented by dark 
matter particles. The aim of our method here is to select a sample 



6 A.R Cooper et al. 



of representative particles from the parent N-body simulation to 
trace each such stellar population, individually. We describe first 
the general objective of our selection process, and then examine 
the selection criteria that we apply in practice. 

Consider first the case of a single galaxy evolving in isola- 
tion. At a given simulation snapshot (B) the total mass of new stars 
formed since the previous snapshot (A) is given by the difference in 
the stellar mass of the semi-analytic galaxy recorded at each time, 

AM^^ = Mf-M^. (1) 

In our terminology, AM^^ is a single stellar population (we do 
not track the small amount of mass lost during subsequent stellar 
evolution). The total mass in metals within the population is deter- 
mined in the same way as the stellar mass; we do not follow in- 
dividual chemical elements. In a similar manner, the luminosity of 
the new population (at 2; = 0) is given by the difference of the total 
luminosities (after evolution to 2; = 0) at successive snapshots. 

From the list of particles in the simulation identified with the 
dark matter halo of the galaxy at B, we select a subset to be trac- 
ers of the stellar population AM^^ . Particles in this tracer set are 
'tagged', i.e. are identified with data describing the stellar popula- 
tion. In the scheme we adopt here, equal 'weight' (fraction of stellar 
mass, luminosity and metals in AM^^) is given to each particle 
in the set of tracers. We repeat this process for all snapshots, apply- 
ing the energy criterion described below to select a new set of DM 
tracers each time new stars are formed in a particular galaxy. In this 
scheme, the same DM particle can be selected as a tracer at more 
than one output time (i.e. the same DM particle can be tagged with 
more than one stellar population). Hence a given DM particle ac- 
cumulates its own individual star formation history. The dynamical 
evolution of satellite haloes determines whether or not a particular 
particle is eligible for the assignment of new stars during any given 
episode of star formation. 

So far we have considered an 'isolated' galaxy. In practice, 
we apply this technique to a merger tree, in which a galaxy grows 
by the accretion of satellites as well as by in situ star formation. 
In the expression above, the total stellar mass at A, M^, is simply 
modified to include a sum over N immediate progenitor galaxies 
in addition to the galaxy itself i.e., 

AM^^ = Mf - M^^o - Mf^, (2) 

i>0 

where M^q represents the galaxy itself and is the total stellar 
mass (at A) of the i'th progenitor deemed to have merged with the 
galaxy in the interval AB. Stars forming in the progenitors during 
the interval AB and stars forming in the galaxy itself are treated as 
a single population. 

There is a one-to-one correspondence between a galaxy and 
a dark matter structure (halo or subhalo) from which particles are 
chosen as tracers of its newly formed stars. As discussed in Sec- 
tion [2]2j a satellite galaxy whose host subhalo is no longer identi- 
fied by SUBFIND has its cold gas content transferred immediately to 
the central galaxy of their common parent halo and forms no new 
stars. In the semi-analytic model, the stars of the satellite are also 
added to the bulge component of the central galaxy. This choice is 
irrelevant in our particle representation, as we can track the actual 
fate of these stars. 



3.2 Assignment criteria 

3. 2. 1 Selection of dark matter particles 

In this section we describe how we choose the dark matter parti- 
cles within haloes that are to be tagged with a newly formed stellar 
population. In Section [T] we briefly described the particle-tagging 
method employed by Bullock & Johnston (2005), the philosophy of 
which we term 'in vitro', using idealised initial conditions to simu- 
late accretion events individually in a 'controlled' environment. By 
contrast, our approach is to postprocess fully cosmological simu- 
lations 'in v/voQ In a fully cosmological N-body simulation the 
growth of the central potential, the structure of the halo and the 
orbits, accretion times and tidal disruption of subhaloes are fully 
consistent with one another. The central potential is non- spherical 
(although no disc component is included in our dynamical model) 
and can grow violently as well as through smooth accretion. Our 
model is therefore applicable at high redshift when the halo is un- 
dergoing rapid assembly. The complexities in the halo potential re- 
alised in a fully cosmological simulation are likely to be an impor- 
tant influence on the dynamics of satellites (e.g. Sales et al. 2007a) 
and on the evolution of streams, through phase-mixing and orbital 
precession (e.g. Helmi & White 1999). 

We approach the selection of dark matter particles for stellar 
tagging differently to Bullock & Johnston (2005), because we are 
postprocessing a cosmological N-body simulation rather than con- 
structing idealised initial conditions for each satellite. Rather than 
assigning the mass-to-light ratio of each tagged particle by compar- 
ing stellar and dark matter energy distribution functions in the halo 
concerned, we assume that the energy distribution of newly formed 
stars traces that of the dark matter. We order the particles in the 
halo by binding energ}!^ and select a most-bound fraction /mb to 
be tagged with newly-formed stars. As previously described, stars 
are shared equally among the selected DM particles. 

Our approach implies a rather simple dynamical model for 
stars in satellite galaxies. However, the main results of this paper 
do not concern the satellites themselves; instead we focus on the 
debris of objects that are totally (or largely) disrupted to build the 
stellar halo. As we describe below, we compare the structure and 
kinematics of our model satellites (those that survive at 2; = 0) to 
Local Group dwarf galaxies in order to fix the value of the free pa- 
rameter, /mb . Since we impose this constraint, our method cannot 
predict these satellite properties ab initio. Constraining our model 
in this way ensures reasonable structural properties in the popula- 
tion of progenitor satellites, and retains full predictive power with 
regard to the stellar halo. More complex models would, of course, 
be possible, in which /mb is not a free parameter but is instead 
physically determined by the semi-analytic model. It would also 
be possible to use a more complicated tagging scheme to attempt 
to represent, for example, star formation in a disc. However, such 
models would add substantial complexity to the method and there 
are currently very few observational constraints on how stars were 
formed in satellite galaxies. Thus, we believe that a simple model 
suffices for our present study of the stellar halo. 

Our approach has similarities with that of De Lucia & Helmi 

^ This terminology should not be taken to imply that 'star particles' them- 
selves are included in the N-body simulation; here stellar populations are 
simply tags affixed to dark matter particles. 

^ Here, the most bound particle is that with the most negative total energy, 
including both kinetic and gravitational contributions. Binding energies are 
computed relative to the bound set of particles comprising an object identi- 
fied by SUBFIND. 
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Figure 2. Examples of individual satellites in our models (solid black lines), compared to Fornax (red) and Carina (blue), showing surface brightness (left, 
Irwin & Hatzidimitriou 1995) and line-of-sight velocity dispersion (right. Walker et al. 2009). With our fiducial GALFORM model, simultaneous matches to 
both (t{R) and iJi{R) for these datasets are found only among satellites that have undergone substantial tidal stripping (see text). 



(2008), who tag the most bound 10% of particles in satellite haloes 
with stars. However, De Lucia & Helmi perform this tagging only 
once for each satellite, at the time at which its parent halo becomes 
a subhalo of the main halo (which we refer to as the time of infalQ). 
Both this approach and that of Bullock & Johnston (2005) define 
the end result of the previous dynamical evolution of an infalling 
satellite, the former by assuming light traces dark matter and the 
latter with a parameterized King profile. 

As described above, in our model each newly-formed stellar 
population is assigned to a subset of DM particles, chosen accord- 
ing to the 'instantaneous' dynamical state of its host halo. This 
choice is independent of any previous star formation in the same 
halo. It is the dynamical evolution of these many tracer sets in each 
satellite that determines its stellar distribution at any point in the 
simulation. 

Implementing a particle-tagging scheme such as this within a 
fully cosmological simulation requires a number of additional is- 
sues to be addressed, which we summarise here. 

(i) Subhalo assignments: Star formation in a satellite galaxy will 
continue to tag particles regardless of the level of its halo in the hi- 
erarchy of bound structures (halo, subhalo, subsubhalo etc.). The 
growth of a dark matter halo ends when it becomes a subhalo of a 
more massive object, whereupon its mass is reduced through tidal 
stripping. The assignment of stars to particles in the central regions 
according to binding energy should, of course, be insensitive to the 
stripping of dark matter at larger radii. However, choosing a fixed 
fraction of dark matter tracer particles to represent new stellar pop- 
ulations couples the mass of the subhalo to the number of particles 
chosen. Therefore, when assigning stars to particles in a subhalo, 
we instead select a fixed number of particles, equal to the number 
constituting the most-bound fraction /mb of the halo at the time of 
infall. 

(ii) Equilibrium criterion: To guard against assigning stars to 
sets of tracer particles that are temporarily far from dynamical 
equilibrium, we adopt the conservative measure of deferring as- 
signments to any halo in which the centres of mass and potential 
are separated by more than 7% of the half-mass radius ri/2. We 

4 In both Bullock & Johnston (2005) and De Lucia & Helmi (2008) only 
satellites directly accreted by the main halo 'trigger' assignments to dark 
matter; the hierarchy of mergers/accretions forming a directly-infalling 
satellite are subsumed in that single assignment. 



select 0.07 ri/2 in accordance with the criterion of 0.14 rvir used 
to select relaxed objects in the study of Neto et al. (2007), taking 
Tvir ^ 2ri/2. These deferred assignments are carried out at the 
next snapshot at which this criterion is satisfied, or at the time of 
infall into a more massive halo. 



(iii) No in situ star formation: Stars formed in the main galaxy 
in each Aquarius simulation (identified as the central galaxy of the 
most massive dark halo at 2; = 0) are never assigned to DM parti- 
cles. This exclusion is applied over the entire history of that galaxy. 
Stars formed in situ are likely to contribute to the innermost regions 
of the stellar halo, within which they may be redistributed in merg- 
ers. However, the dynamics of stars formed in a dissipationally- 
collapsed, baryon-dominated thin disc cannot be represented with 
particles chosen from a dark matter-only simulation. We choose in- 
stead to study the accreted component in isolation. Our technique 
none the less offers the possibility of extracting some information 
on a fraction of in situ stars were we to assign them to dark matter 
particles (those contributing to the bulge or forming at early times, 
for example). We choose to omit this additional complexity here. 
SPH simulations of stellar haloes (which naturally model the in 
situ component more accurately than the accreted component) sug- 
gest that the contribution of in situ stars to the halo is small beyond 
- 20 kpc (Abadi et al. 2006; Zolotov et al. 2009). 

At early times, when the principal halo in each simulation is 
growing rapidly and near-equal-mass mergers are common, the def- 
inition of the 'main' branch of its merger tree can become ambigu- 
ous. Also, the main branch of the galaxy merger tree need not fol- 
low the main branch of the halo tree. Hence, our choice of which 
branch to exclude (on the basis that it is forming 'in situ' stars) also 
becomes ambiguous; indeed, it is not clear that any of these 'equiv- 
alent' early branches should be excluded. Later we will show that 
two of our haloes have concentrated density profiles. We have con- 
firmed that these do not arise from making the 'wrong' choice in 
these uncertain cases, i.e. from tagging particles in the dynamically 
robust core of the 'true' main halo. Making a different choice of 
the excluded branch in these cases (before the principal branch can 
be unambiguously identified) simply replaces one of these concen- 
trated components with another very similar component. Therefore, 
we adopt the above definition of the galaxy main branch when ex- 
cluding in situ stars. 
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3.2.2 Individual satellites 

We show in the following section that with a suitable choice of the 
most-bound fraction, our method produces a population of model 
satellites at 2; = having properties consistent with observed re- 
lationships between magnitude, half-light radius/surface brightness 
and velocity dispersion for satellite populations of the Milky Way 
and M31. In Fig. [2 we show profiles of surface brightness and ve- 
locity dispersion for two individual satellites from these models 
at ^ = 0, chosen to give a rough match to observations of For- 
nax and Carina. This suggests that our galaxy formation model and 
the simple prescription for the spatial distribution of star formation 
can produce realistic stellar structures within dark haloes. How- 
ever, while it is possible to match these individual observed satel- 
lites with examples drawn from our models, we caution that we can 
only match their observed surface brightness and velocity disper- 
sion profiles simultaneously by choosing model satellites that have 
suffered substantial tidal stripping. This is most notable in the case 
of our match to Fornax, which retains only 2% of its dark matter 
relative to the time of its accretion to the main halo, and 20% of 
its stellar mass. However, as we show in Section I4I21 the majority 
of massive surviving satellites have not suffered substantial tidal 
stripping. 

We have tested our method with assignments for each satellite 
delayed until the time of infall, as in De Lucia & Helmi (2008). 
This results in slightly more compact galaxies than in our standard 
in vivo approach, where mergers and tidal forces (and relaxation 
through two-body encounters for objects near the resolution limit) 
can increase the energies of tagged dark matter particles. However, 
we find that this makes little difference to the results that we discuss 
below. 



3.2.3 Parameter constraints and convergence 

We now compare the 2; = satellite populations of our models 
with trends observed in the dwarf companions of the Milky Way 
and M3 1 in order to determine a suitable choice of the fixed frac- 
tion, /mb, of the most bound dark matter particles selected in a 
given halo. Our aim is to study the stellar halo, and therefore we 
use the sizes of our surviving satellites as a constraint on /mb and 
as a test of convergence. Within the range of /mb that produces 
plausible satellites, the gross properties of our haloes, such as total 
luminosity, change by only a few percent. 

In Fig.O we show the relationship between the absolute mag- 
nitudes, Mv, of satellites (combining data from two of our simu- 
lations, Aq-A and Aq-F), and the projected radius enclosing one 
half of their total luminosity, which we refer to as the effective ra- 
dius, Teff . We compare our models to a compilation of dwarf galaxy 
data in the Local Group, including the satellites of the Milky Way 
and M3 1 . The slope of the median relation for our satellites agrees 
well with that of the data for the choices /mb = 1% and 3%. It 
is clear that a choice of 5% produces bright satellites that are too 
extended, while for 0.5% they are too compact. We therefore prefer 
/mb = 1%. A more detailed comparison to the data at this level 
is problematic: the observed sample of dwarf galaxies available at 
any given magnitude is small, and the data themselves contain puz- 
zling features such as an apparently systematic difference in size 
between the bright Milky Way and M31 satellites. 

Fig. [3] also shows (as dotted lines) the same results for our 
model run on the lower-resolution simulations of haloes Aq-A and 
Aq-F. The particle mass in the Aq-3 series is approximately three 
times greater than in Aq-2, and the force softening scale is larger 
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Figure 3. Median effective radius reff (enclosing half of the total luminos- 
ity in projection) as a function of magnitude for model satellites in haloes 
Aq-A and Aq-F at 2; = 0. A thin vertical dashed line indicates the softening 
scale of the simulation: rgff is unreliable close to this value and meaning- 
less below it. Thick lines represent our higher-resolution simulations (Aq-2) 
using a range of values of the fraction of most bound particles chosen in a 
stellar population assignment, /mb • Dotted lines correspond to lower res- 
olution simulations (Aq-3) of the same haloes. A thick dashed line shows 
the corresponding median of observations of Local Group dwarf galaxies. 
These galaxies, and our model data points for all haloes in the Aq-2 series 
with /mb = 1%, are plotted individually in Fig.|4] 



by a factor of two. We concentrate on the convergence behaviour 
of our simulations for galaxies larger than the softening length, and 
also where our sample provides a statistically meaningful number 
of galaxies at a given magnitude; this selection corresponds closely 
to the regime of the brighter dwarf spheroidal satellites of the Milky 
Way and M31, —15 < My < —5. In this regime. Fig. [3] shows 
convergence of the median relations brighter than My = —5 for 
/mb = 3% and 5%. The case for /mb = 1% is less clear-cut. 
The number of particles available for a given assignment is set by 
the mass of the halo; haloes near the resolution limit (with 100 
particles) will, of course, have only ^ 1 particle selected in a sin- 
gle assignment. In addition to this poor resolution, galaxies formed 
by such small-number assignments are more sensitive to spurious 
two-body heating in the innermost regions of subhaloes. We there- 
fore expect the resulting galaxies to be dominated by few-particle 
'noise' and to show poor convergence behaviour. 

We adopt /mb = 1% as a reasonable match to the data (noting 
also that it lies close to the power-law fit employed by Bullock & 
Johnston (2005) to map luminosities to satellite sizes). We believe 
the resulting satellites to be sufficiently converged at the resolution 
of our Aq-2 simulations with this choice of /mb to permit a sta- 
tistical study of the disrupted population represented by the stellar 
halo. In support of this assertion, we offer the following heuris- 
tic argument. The change in resolution from Aq-3 to Aq-2 results 
in approximately three times more particles being selected at fixed 
/mb; likewise, a change in /mb from 1% to 3% selects three times 
more particles at fixed resolution. Therefore, as /mb = 3% has 
converged at the resolution of Aq-3, it is reasonable to expect that 
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Figure 4. Projected half-light radius (left), mean luminosity-weighted ID velocity dispersion (centre) and central surface brightness (right) of simulated 
satellite galaxies (defined by roc < 280 kpc) that survive in all haloes at 2; = 0, as a function of absolute V-band magnitude. Observational data for Milky 
Way and M31 satellites are shown as orange symbols; values are from Mateo (1998) and other authors as follows: bright satellites (triangles pointing right, 
Grebel, Gallagher & Harbeck 2003); faint MW satellites discovered since 2005 (triangles pointing up, Martin, de Jong & Rix 2008); M31 dwarf spheroidals 
(triangles pointing left, McConnachie et al. 2006; Martin et al. 2009); M31 ellipticals (squares); Local Group 'field' dwarf spheroidals and dwarf irregulars 
(stars). In the central panel we use data for Milky Way satellites only tabulated by Wolf et al. (2009) and for the SMC, Grebel et al. (2003). In the rightmost 
panel, we plot data for the Milky Way and M31 (Grebel et al. 2003; Martin et al. 2008). A dashed line indicates the surface brightness of an object of a given 
magnitude with reff = 2.8e, the gravitational softening scale (see Section [Zll . 



/mb = 1% selects a sufficient number of particles to ensure that 
satellite sizes are not dominated by noise at the resolution of Aq-2. 
We show below that the most significant contribution to the halo 
comes from a handful of well resolved objects with My < —10, 
rather than from the aggregation of many fainter satellites. Addi- 
tionally, as demonstrated for example by Penarrubia, McConnachie 
& Navarro (2008a); Penarrubia, Navarro, & McConnachie (2008b); 
Penarrubia et al. (2009), there is a 'knife-edge' between the onset 
of stellar stripping and total disruption for stars deeply embedded 
within the innermost few percent of the dark matter in a halo. We 
conclude that premature stripping resulting from an over-extension 
of very small satellites in our model is unlikely to alter the gross 
properties of our stellar haloes. 

The points raised above in connection with Fig. [3] make clear 
that the in vivo particle tagging approach demands extremely high 
resolution, near the limits of current cosmological N-body simula- 
tions. The choice of /mb = 1% in this approach (from an accept- 
able range of 1 — 3%) is not arbitrary. For example, a choice of 
/mb = 10% (either as a round-number estimate. 

For the remainder of this paper we concentrate on the higher 
resolution Aq-2 simulations. In Fig.|4lwe fix /mb at 1% and com- 
pare the surviving satellites from all six of our haloes with observa- 
tional data for three properties correlated with absolute magnitude: 
effective radius, reff, mean luminosity- weighted line-of-sight ve- 
locity dispersion, a, and central surface brightness, /xo (although 
the latter is not independent of reff). In all cases our model satel- 
lites agree well with the trends and scatter in the data brighter than 
Mv -5. 

The force softening scale of the simulation (indicated in the 
first and third panels by dashed lines) effectively imposes a maxi- 
mum density on satellite dark haloes. At this radial scale we would 
expect reff to become independent of magnitude for numerical rea- 
sons: Fig. |4] shows that the reff (My) relation becomes steeper for 
galaxies fainter than My ^ —9 , corresponding to reff ^ 200 pc. 



This resolution-dependent maximum density corresponds to a min- 
imum surface brightness at a given magnitude. The low-surface- 
brightness limit in the Milky Way data shown in the right-hand 
panel of Fig. |4l corresponds to the completeness limit of current 
surveys (e.g. Koposov et al. 2008; Tollerud et al. 2008). The lower 
surface brightness satellite population predicted by our model is 
not, in principle, incompatible with current data. 

In Fig. [5] we show the relationship between total luminosity 
and the mass of dark matter enclosed within 300 pc, M300, for our 
simulated satellites in all haloes. This radial scale is well-resolved 
in the level 2 Aquarius simulations (see also Font et al. 2009, in 
prep.). Our galaxies show a steeper trend than the data of Strigari 
et al. (2008), with the strongest discrepency (0.5 dex in M300) for 
the brightest satellites. Nevertheless, both show very little variation, 
having M300 10^ Mq over five orders of magnitude in luminos- 
ity. In agreement with previous studies using semi-analytic models 
and lower-resolution N-body simulations (Maccio, Kang & Moore 
2009; Busha et al. 2009; Li et al. 2009b; Koposov et al. 2009), and 
N-body gasdynamic simulations (Okamoto & Frenk 2009), we find 
that this characteristic scale arises naturally as a result of astrophys- 
ical processes including gas cooling, star formation and feedback. 

3.3 Defining the stellar halo and satellite galaxies 

To conclude this section, we summarise the terminology we adopt 
when describing our results. Tagged dark matter particles in the 
self-bound haloes and subhaloes identified by SUB FIND consti- 
tute our 'galaxies'. Our stellar haloes comprise all tagged parti- 
cles bound to the main halo in the simulation, along with those 
tagged particles not in any bound group (below we impose an ad- 
ditional radial criterion on our definition of the stellar halo). All 
galaxies within 280 kpc of the centre of the main halo are classed 
as 'satellites', as in the luminosity functions shown in Fig.[T] Cen- 
tres of mass of the stellar haloes and satellites are determined from 
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Figure 5. Mass in dark matter enclosed within 300 pc (M300) as a func- 
tion of luminosity (V-band) for satellites in each of our simulated haloes 
(coloured points, colours as Fig.[TJ. Maximum likelihood values of M300 
for Milky Way dwarf spheroidals from Strigari et al. (2008) are shown (or- 
ange squares), with error bars indicating the range with likelihood greater 
than 60.6% of the maximum. 



tagged particles only, using the iterative centring process described 
by Power et al. (2003). 

Many structural elements of a galaxy intermix within a few 
kiloparsecs of its centre, and attempts to describe the innermost 
regions of a stellar halo require a careful and unambiguous defini- 
tion of other components present. This is especially important when 
distinguishing between those components that are represented in 
our model and those that are not. Therefore, before describing our 
haloefl we first summarise some of these possible sources of con- 
fusion, clarify what is and is not included in our model, and define a 
range of galactocentric distances on which we will focus our anal- 
ysis of the stellar halo. 

As discussed above, our model does not track with parti- 
cles any stars formed in situ in the central 'Milky Way' galaxy, 
whether in a rotationally supported thin disc or otherwise (this cen- 
tral galaxy is, of course, included in the underlying semi-analytic 
model). We therefore refer to the halo stars that are included in our 
model as accreted and those that form in the central galaxy (and 
hence are not explicitly tracked in our model) as in situ. Observa- 
tional definitions of the 'stellar halo' typically do not attempt to 
distinguish between accreted and in situ stars, only between com- 
ponents separated empirically by their kinematic, spatial and chem- 
ical properties. 

The 'contamination' of a purely- accreted halo by stars formed 
in situ is likely to be most acute near the plane of the disc. Ob- 
servations of the Milky Way and analogous galaxies frequently 
distinguish a 'thick disc' component (Gilmore & Reid 1983; Car- 
ollo et al. 2009) thought to result either from dynamical heating of 
the thin disc by minor mergers (e.g. Toth & Ostriker 1992; Quinn, 

^ We explicitly distinguish between the stellar halo and the dark halo in 
ambiguous cases; typically the former is implied throughout 



Hernquist & Fullagar 1993; Velazquez & White 1999; Font et al. 
2001; Benson et al. 2004; Kazantzidis et al. 2008) or from accre- 
tion debris (Abadi et al. 2003; Yoachim & Dalcanton 2005, 2008). 
The presence of such a component in M31 is unclear: an 'extended 
disc' is observed (Ibata et al. 2005), which rotates rapidly, con- 
tains a young stellar population and is aligned with the axes of the 
thin disc, but extends to 40 kpc and shows many irregular mor- 
phological features suggestive of a violent origin. In principle, our 
model will follow the formation of accreted thick discs. However, 
the stars in our model only feel the potential of the dark halo; the 
presence of a massive baryonic disc could significantly alter this 
potential in the central region and influence the formation of an 
accreted thick disc (e.g. Velazquez & White 1999). 

Our models include that part of the galactic bulge built from 
accreted stars, but none of the many other possible processes of 
bulge formation (starbursts, bars etc.). However, the interpretation 
of this component, the signatures of an observational counterpart 
and the extent to which our simulation accurately represents its dy- 
namics are all beyond the scope of this paper. Instead, we will con- 
sider stars within 3 kpc of the dark halo potential centre as 'accreted 
bulge', and define those between 3 kpc and a maximum radius of 
280 kpc as the 'stellar halo' on which we will focus our analysis. 
This arbitrary radial cut is chosen to exclude the region in which the 
observational separation of 'bulge' and 'halo' stars is not straight- 
forward, and which is implicitly excluded from conventional ob- 
servational definitions of the halo. It is not intended to reflect a 
physical scale-length or dichotomy in our stellar haloes, analogous 
to that claimed for the Milky Way (e.g. Carollo et al. 2007, 2009). 
Beyond 3 kpc we believe that the ambiguities discussed above and 
the 'incompleteness' of our models with regard to stars formed in 
situ should not substantially affect the comparison of our accreted 
stars with observational data. 



4 RESULTS: THE AQUARIUS STELLAR HALOES 

In this section, we present the six stellar haloes resulting from the 
application of the method described above to the Aquarius simu- 
lations. Here our aim is to characterise the assembly history of the 
six haloes and their global properties. Quantities measured for each 
halo are collected in Tabled These include a measure of the num- 
ber of progenitor galaxies contributing to the stellar halo, A^prog. 
This last quantity is not the total number of accreted satellites, but 
instead is defined as A^prog = ^haio/ ^prog,i where mprog,i is 
the stellar mass contributed by the i'th progenitor. A^prog is equal to 
the total number of progenitors in the case where each contributes 
equal mass, or to the number of significant progenitors in the case 
where the remainder provide a negligible contribution. 



4.1 Visualisation in projection 

A 300 X 300 kpc projected surface brightness map of each stellar 
halo at z = is shown in Fig. [6] Substantial diversity among the 
six haloes is apparent. Haloes Aq-B and Aq-E are distinguished by 
their strong central concentration, with few features of detectable 
surface brightness beyond ^ 20 kpc. Haloes Aq-A, Aq-C, Aq-D 
and Aq-F all show more extended envelopes to 75-100 kpc; each 
envelope is a superposition of streams and shells that have been 
phase-mixed to varying degrees. 

Analogues of many morphological features observed in the 
halo of M31 (Ibata et al. 2007; Tanaka et al. 2009; McConnachie 
et al. 2009) and other galaxies (e.g. Martmez-Delgado et al. 2008) 
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Figure 6. V-band surface brightness of our model haloes (and surviving satellites), to a limiting depth of 35 mag/arcsec^. The axis scales are in kiloparsecs. 
Only stars formed in satellites are present in our particle model; there is no contribution to these maps from a central galactic disc or bulge formed in situ (see 
Section (33) 
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Figure 7. As Fig.|6] but here showing only those stars stripped from satellites that survive at z = 
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Table 2. For each of our simulated haloes we tabulate: the luminosity and mass of halo stars (in the range 3 < r < 280 kpc); the mass of 
accreted bulge stars (r < 3 kpc); the total stellar mass and V-band magnitude of the central galaxy in GALFORM; the number of surviving 
satellites (brighter than My = 0); the fraction of the total stellar mass within 280 kpc bound in surviving satellites at 2; = 0, /sat; the fraction 
of halo stellar mass (r < 280 kpc) contributed by these surviving satellites, /surv; the number of halo progenitors, A^prog (see text); the 
half-light radius of the stellar halo (r < 280 kpc); the inner and outer slope and break radius of a broken power-law fit to the three-dimensional 
density profile of halo stars (3 < r < 280 kpc). 



Halo ^y,halo ^★,halo ^★jbulge ^gal ^sat /sat /surv A^prog ^1/2 n '^out ^brk 

[IO^Lq] [IO^Mq] [IO^Mq] [IOIOMq] [kpc] [kpc] 



A 


1.51 


2.80 


1.00 


1.88 


-20.3 


161 


0.61 


0.065 


3.8 


20 


-2.7 


-8.2 


80.4 


B 


1.27 


2.27 


3.33 


1.49 


-20.1 


91 


0.07 


0.036 


2.4 


2.3 


-4.2 


-5.8 


34.6 


C 


1.95 


3.58 


0.34 


7.84 


-21.3 


150 


0.28 


0.667 


2.8 


53 


-2.0 


-9.4 


90.8 


D 


5.55 


9.81 


1.32 


0.72 


-19.1 


178 


0.35 


0.620 


4.3 


26 


-2.0 


-5.9 


37.7 


E 


0.90 


1.76 


16.80 


0.45 


-18.6 


135 


0.11 


0.003 


1.2 


1.0 


-4.7 


-4.4 


15.2 


F 


17.34 


24.90 


6.42 


1.36 


-20.1 


134 


0.28 


0.002 


1.1 


6.3 


-2.9 


-5.9 


14.0 



can be found in our simulations. For example, the lower left quad- 
rant of Aq-A shows arc-like features reminiscent of a complex 
of 'parallel' streams in the M31 halo labelled A, B, C and D by 
Ibata et al. (2007) and Chapman et al. (2008), which have surface 
brightnesses of 30 — 33 mag arcsec"^ and a range of metallici- 
ties (Tanaka et al. 2009). These streams in Aq-A can also be traced 
faintly in the upper right quadrant of the image and superficially re- 
semble the edges of 'shells'. In fact, they result from two separate 
progenitor streams, each tracing multiple wraps of decaying orbits 
(and hence contributing more than one 'arc' each). Seen in three di- 
mensions, these two debris complexes (which are among the most 
significant contributors to the Aq-A halo) are elaborate and irreg- 
ular structures, the true nature of which is not readily apparent in 
any given projectior|fl 

The brightest and most coherent structures visible in Fig. [6l 
are attributable to the most recent accretion events. To illustrate 
the contribution of recently-infalling objects (quantified in Sec- 
tion |4.21 ), we show the same projections of the haloes in Fig.[7l but 
include only those stars whose parent satellite survives at z = 0. 
In haloes Aq-C and Aq-D, stars stripped from surviving satellites 
constitute 60 — 70% of the halo, while in the other haloes their 
contribution is < 10%. Not all the recently-infalling satellites re- 
sponsible for bright halo features survive; for example, the massive 
satellite that merges at z ^ 0.3 and produces the prominent set of 
'shells' in Aq-F. 

Fig. [6l shows that all our haloes are notably flattened, particu- 
larly in the central regions where most of their light is concentrated. 
Axial ratios q — c/a and s — b/aof three-dimensional ellipsoidal 
fits to halo stars within 10 kpc of the halo centre are given in Table[3] 
(these fits include stars within the accreted bulge region defined 
above). Most of our haloes are strongly prolate within 10 kpc. Halo 
Aq-E is very different, having a highly oblate (i.e. disc-like) shape 
in this region - this structure of ^ 20 kpc extent can be seen 'edge 
on' in Fig.[6]and can be described as an 'accreted thick disc' (e.g. 
Abadi et al. 2003; Penarrubia, McConnachie & Babul 2006; Read 
et al. 2008). We defer further analysis of this interesting object to a 
subsequent paper. Beyond 10-30 kpc, the stellar mass in our haloes 
is not smoothly distributed but instead consists of a number of dis- 
crete streams, plumes and other irregular structures. Fits to all halo 
stars assuming a smoothly varying ellipsoidal distribution of mass 



° Three orthogonal projections for each halo can be found at 
|http : / / www . virgo . dur . ac . uk/ aquarius] 



Table 3. Axial ratios q = c/a and s = 6/a of stellar-mass-weighted three- 
dimensional ellipsoidal fits to halo stars within a galactocentric radius of 
10 kpc. These were determined using the iterative procedure described by 
AUgood et al. (2006), which attempts to fit the shapes of self-consistent 
'isodensity' contours. A spherical contour of r = 10 kpc is assumed ini- 
tially; the shape and orientation of this contour are then updated on each 
iteration to those obtained by diagonialzing the inertia tensor of the mass 
enclosed (maintaining the length of the longest axis). The values thus ob- 
tained are slightly more prolate than those obtained from a single diagnon- 
alization using all mass with a spherical contour (i.e. the first iteration of 
our approach), reflecting the extremely flattened shapes of our haloes at this 
radius. The oblate shape of Aq-E is not sensitive to this choice of method. 



Halo 


A 


B 


C 


D 


E 


F 


qio 
sio 


0.27 
0.30 


0.28 
0.32 


0.29 
0.32 


0.33 
0.42 


0.36 
0.96 


0.21 
0.25 



interior to a given radius do not accurately describe these sparse 
outer regions. 

Few observations of stellar halo shapes are available for com- 
parison with our models. M31 is the only galaxy in which a pro- 
jected stellar halo has been imaged to a depth sufficient to ac- 
count for a significant fraction of halo stars. Pritchet & van den 
Bergh (1994) measured a projected axial ratio for the M31 halo at 
^ 10 kpc of ^ 0.5. Ibata et al. (2005) describe a highly irregular 
and rotating inner halo component or 'extended disc' (to ^ 40 kpc) 
of 27 — 31 mag/arcsec^, aligned with the thin disc and having an 
axial ratio ^ 0.6 in projection. Zibetti & Ferguson (2004) find a 
similar axial ratio for the halo of a galaxy at z = 0.32 observed 
in the Hubble ultra-deep field. Evidence for the universality of flat- 
tened stellar haloes is given by Zibetti, White & Brinkmann (2004), 
who find a best-fitting projected axial ratio of 0.5 — 0.7 for the 
low surface brightness envelope of ^ 1000 stacked edge-on late- 
type galaxies in SDSS. A mildly oblate halo with c/a ^ 0.6 is re- 
ported for the Milky Way, with an increase in flattening at smaller 
radii « 20 kpc; e.g. Chiba & Beers 2000; Bell et al. 2008; Car- 
olio et al. 2007). Interestingly, Morrison et al. (2009) present evi- 
dence for a highly flattened halo (c/a ^ 0.2) component in the So- 
lar neighbourhood, which appears to be dispersion- supported (i.e. 
kinematically distinct from a rotationally supported thick disc). 

The shapes of components in our haloes selected by their kine- 
matics, chemistry or photometry may be very different to those 
obtained from the aggregated stellar mass. A full comparison, ac- 
counting for the variety of observational selections, projection ef- 
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Figure 8. The growth of the stellar halo (upper panel) and the dark matter 
halo (the principal branch; lower panel) as a function of expansion factor 
(bottom axis) or redshift (top axis). Lines show the mass fraction of each 
halo in place at a given time. Stars are counted as belonging to the stellar 
halo when the DM particle that they tag is assigned to the principal halo, or 
is not bound to any SUB find group. 



fects and definitions of 'shape' used in the measurements cited 
above, is beyond the scope of this paper. We emphasize, however, 
that the flattening in our stellar haloes cannot be attributed to any 
'baryonic' effects such as a thin disc potential (e.g. Chiba & Beers 
2001) or star formation in dissipative mergers and bulk gas flows 
(e.g. Bekki & Chiba 2001). Furthermore, it is unlikely to be the re- 
sult of a (lesser) degree of flattening in the dark halo. Instead the 
structure of these components is most likely to reflect the intrin- 
sically anisotropic distribution of satellite orbits. In certain cases 
(for example, Aq-D and Aq-A), it is clear that several contributing 
satellites with correlated trajectories are responsible for reinforcing 
the flattening of the inner halo. 



4.2 Assembly history of the stellar halo 

We now examine when and how our stellar haloes were assembled. 
Fig. [8] shows the mass fraction of each stellar halo (here including 
the accreted bulge component defined in Section (33]) in place (i.e. 
unbound from its parent galaxy) at a given redshift. We count as 
belonging to the stellar halo all 'star particles' bound to the main 
dark halo and within 280 kpc of its centre at 2; = 0. This is com- 
pared with the growth of the corresponding host dark haloes. Our 
sample spans a range of assembly histories for haloes even though 
the halos have very similar final mass. 

Not surprisingly, the growth of the dark halo is considerably 
more smooth than that of the stellar halo. The 'luminous' satellite 
accretion events contributing stars are a small subset of those that 
contribute to the dark halo, which additionally accretes a substantial 
fraction of its mass in the form of 'diffuse' dark matter (Wang et al. 
in prep.). As described in detail by Penarrubia et al. (2008a,b), the 
dark haloes of infalling satellites must be heavily stripped before 
the deeply embedded stars are removed. This gives rise to time-lags 
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Figure 9. Luminosity functions of surviving satellites (solid) in each of 
our six haloes, compared with those of totally disrupted halo progenitors 
(dashed). These are constructed using only stars formed in each satellite 
before the time of infall (the halo-subhalo transition). The luminosity of 
each population is that after evolution to 2; = 0. 



seen in Fig. [8] between the major events building dark and stellar 
haloes. 

To characterise the similarities and differences between their 
histories, we subdivide our sample of six stellar haloes into two 
broad categories: those that grow through the gradual accretion of 
many progenitors (Aq-A, Aq-C and Aq-D) and those for which 
the majority of stellar mass is contributed by only one or two ma- 
jor events (Aq-B, Aq-E and Aq-F). We refer to this latter case as 
'few-progenitor' growth. The measure of the number of 'most- 
significant' progenitors given in Table [2l A^prog, also ranks the 
haloes by the 'smoothness' of their accretion history, reflecting the 
intrinsically stochastic nature of their assembly. 

Fig. [9] compares the luminosity functions (LFs) of surviving 
satellites with that of those totally disrupted to form the stellar halo, 
measuring luminosity at the time of infall in both cases. In general, 
there are fewer disrupted satellites than survivors over almost all 
luminosities, although the numbers and luminosities of the very 
brightest contributors and survivors are comparable in each halo. 
The deficit in the number of disrupted satellites relative to survivors 
is most pronounced in the few-progenitor haloes Aq-B and Aq-F. 

Fig. [TOl summarises the individual accretion events that con- 
tribute to the assembly of the stellar halo, plotting the stellar mass 
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Figure 10. Main panel: for satellites that have been stripped to form the 
stellar haloes, symbols show the redshift of infall and total mass contributed 
to the stellar halo at z = (in the range 3 < r < 280 kpc). Vertical 
lines indicate the total mass of each stellar halo in this radial range. The 
right-hand y-axis is labelled by lookback time in gigayears. We plot only 
those satellites whose individual contributions, accumulated in rank order 
from the most significant contributor, account for 95% of the total stellar 
halo mass. Satellites totally disrupted by 2; = are plotted as open circles, 
surviving satellites as filled squares (in almost all cases the contributions of 
these survivors are close to their total stellar masses; see text). Lower panel: 
symbols indicate the approximate masses of bright MW satellites, assuming 
a stellar mass-to-light ratio of 2; the Sgr present-day mass estimate is that 
given by Law, Johnston & Majewski (2005). The shaded region indicates 
an approximate range for the MW halo mass in our halo regime (see e.g. 
Bell et al. 2008). 



of the most significant progenitor satellites against their redshift of 
infall (the time at which their host halo first becomes a subhalo 
of the main FOF group). Here we class as significant those satel- 
lites which together contribute 95% of the total halo stellar mass 
(this total is shown as a vertical line for each halo) when accumu- 
lated in rank order of their contribution. By this measure there are 
(5,6,8,6,6,1) significant progenitors for haloes (A,B,C,D,E,F). We 
also compare the masses of the brightest Milky Way satellites to the 
significant contributors in our stellar haloes. Typically the most sig- 
nificant contributors have masses comparable to the most massive 
surviving dwarf spheroidals, Fornax and Sagittarius. 

With the exception of Aq-F, all the most significant contribu- 
tors to our stellar haloes were accreted more than 8 Gyr ago. We 
highlight (as filled squares) those contributors whose cores survive 
as self-bound objects at z = 0. We find that surviving satellites 
accreted before 2; = 1 are the dominant contributors to the many- 
progenitor haloes Aq-C and Aq-D. The extreme case of Aq-F is 
atypical: more than 95% of the halo was contributed by the late 
merger of an object of stellar mass greater than the SMC infalling at 
z ^ 0.7, which does not survive. By contrast, the two least massive 
haloes Aq-B and Aq-E are built by many less massive accretions at 
higher redshift, with surviving satellites making only a minor con- 
tribution (< 10%). Halo Aq-A represents an intermediate case, in 
which stars stripped from a relatively late-infalling survivor add 
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Figure 11. Cumulative mass fraction of each stellar halo originating in 
satellites of stellar mass less than Mgat- Satellite masses are normalised 
to the total stellar halo mass M^aio in ^^ch case, as defined in Section [33l 



significantly 10%) to the mass of a halo predominantly assem- 
bled at high redshift. The relative contributions to the halo of all 
accretion events are illustrated in Fig. \TT\ Each line in this figure 
indicates the fraction of the total halo stellar mass that was con- 
tributed by satellites donating less than a given fraction of this total 
individually. An interesting feature illustrated by this figure con- 
cerns Aq-B, one of our few-progenitor haloes (shown as light blue 
in all figures). Although Fig. [8] shows that the assembly of this halo 
proceeds over time by a series of concentrated 'jumps' in mass, its 
final composition is even less biased to the most significant progen- 
itor than any of the many-progenitor haloes. 

In general, surviving contributors to the halo retain less than 
5% of the total stellar mass that formed in them. A small number of 
surviving contributors retain a significant fraction of their mass, for 
example, the surviving contributor to Aq-A, which retains 25%. In 
Fig. [T2I we show histograms of the number of all surviving satel- 
lites (combining all six haloes) that have been stripped of a given 
fraction of their mass. Most satellites are either largely unaffected 
or almost totally stripped, indicating that the time spent in an inter- 
mediate disrupting state is relatively short. 

In Tabled we give the fraction of mass in the stellar halo that 
has been stripped from surviving satellites, /surv As previously 
stated, this contribution is dominant in haloes Aq-C (67%) and Aq- 
D (62%), significant in Aq-A (7%) and Aq-B (4%) and negligible 
in Aq-E and Aq-F Sales et al. (2007b) find that only - 6% of 
stars in the eight haloes formed in the SPH simulations of Abadi 
et al. (2006) are associated with a surviving satellite. The lack of 
surviving satellites may be attributable to the limited resolution of 
those simulations; clearly, the number of 'survivors' is sensitive to 
the lowest mass at which remnant cores can be resolved. However, 
Bullock & Johnston (2005), and the companion study of Font et al. 
(2006), also conclude that the contribution of surviving satellites is 
small (< 10% in all of their 11 haloes and typically < 1%). As the 
resolution of their simulations is comparable to ours, the predomi- 
nance of surviving contributors in two of our haloes is significant. 
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Figure 12. Number of surviving satellites (aggregated over all six haloes) 
which have lost a fraction, /stripped, of the stellar mass through tidal strip- 
ping. Satellites are divided into three mass bins: massive (purple), inter- 
mediate (dashed orange) and low-mass (dotted black) as quantified in the 
legend. The leftmost bin (demarcated by a vertical line) shows the number 
of satellites that have not suffered any stellar mass loss. 



Bullock & Johnston find that their haloes are built from a sim- 
ilar (small) number of massive objects to ours (e.g. figure 10 of 
Bullock & Johnston 2005) with comparable accretion times (> 8 
Gyr), suggesting that there are no fundamental differences in the 
infall times and masses of accreted satellites. Notably, Font et al. 
(2006) observe that no satellites accreted > 9 Gyr ago survive in 
their subsample of four of the Bullock & Johnston haloes, whereas 
we find that some satellites infalling even at redshifts z > 2 may 
survive (see also Fig. [16] below). The discrepancy appears to stem 
from the greater resilience of satellites accreted at z > 1 in our 
models, including some which contribute significantly to the stel- 
lar haloes. In other words, our model does not predict any more 
late-infalling contributors than the models of Bullock & Johnston. 
The more rapid disruption of massive subhaloes in the Bullock & 
Johnston models may be attributable to one or both of the analytic 
prescriptions employed by those authors to model the growth of 
the dark matter halo and dynamical friction in the absence of a live 
halo. It is also possible that the relation between halo mass and 
concentration assumed in the Bullock & Johnston model results in 
satellites that are less concentrated than subhaloes in the Aquarius 
simulations. 

Current observational estimates (e.g. Bell et al. 2008) imply 
that the stellar halo of the Milky Way is intermediate in mass be- 
tween our haloes Aq-C and Aq-D; if its accretion history is, in fact, 
qualitatively similar to these many-progenitor haloes. Fig. [Tol im- 
plies that it is likely to have accreted its four or five most significant 
contributors around z ^ 1 — 3 in the form of objects with masses 
similar to the Fornax or Leo I dwarf spheroidals. Between one and 
three of the most recently accreted, and hence most massive con- 
tributors, are expected to retain a surviving core, and to have a stel- 
lar mass comparable to Sagittarius (Msgr 5 x 10^ Mq or ^ 50% 



of the totaQ halo mass, infalling at a lookback time of 5 Gyr; 
Law et al. 2005). It is also possible that the Canis Major overden- 
sity (with a core luminosity comparable to that of Sagittarius; Mar- 
tin et al. 2004) associated with the low-latitude Monoceros stream 
(Newberg et al. 2002; Yanny et al. 2003; Ibata et al. 2003) should 
be included in the census of 'surviving contributors' (although this 
association is by no means certain; e.g. Mateu et al. 2009). There- 
fore, the picture so far established for the Milky Way appears to be 
in qualitative agreement with the presence of surviving cores from 
massive stellar halo contributors in our simulations. 



4.3 Bulk halo properties and observables 

4. 3. 1 Distribution of mass 

In Fig.[T3]we show the spherically averaged density profiles of halo 
stars (excluding material bound in surviving satellites, but making 
no distinction between streams, tidal tails or other overdensities, 
and a 'smooth' component). The notable degree of substructure in 
these profiles contrasts with the smooth dark matter haloes, which 
are well-fit by the Einasto profiles shown in Fig. [13] As discussed 
further below, this stellar substructure is due to the contribution 
of localised, spatially coherent subcomponents within the haloes, 
which are well resolved in our particle representation. 

The shapes of the density profiles are broadly similar, show- 
ing a strong central concentration and an outer decline consider- 
ably steeper than that of the dark matter. We overplot in Fig. [13] 
an approximation of the Milky Way halo profile (Bell et al. 2008) 
and normalization (Fuchs & JahreiB 1998; Gould et al. 1998). The 
gross structure of our three many-progenitor haloes Aq-A, Aq-C 
and Aq-D can be fit with broken power-law profiles having indices 
similar to the Milky Way {n ^ —3) interior to the break. Bell et al. 
(2008) note that their best-fitting observational profiles do not fully 
represent the complex structure of the halo, even though they mask 
out known overdensities (our fits include all halo substructure). Our 
fits decline somewhat more steeply than the Bell et al. data beyond 
their break radii. We suggest that the Milky Way fit may represent 
variation at the level of the fluctuations seen in our profiles, and that 
an even steeper decline may be observed with a representative and 
well- sampled tracer population to > 100 kpc (For example, Ivezic 
et al. 2000, find a sharp decline in counts of RR Lyr stars beyond 
^ 60 kpc). In contrast with the many-progenitor haloes, two of our 
few-progenitor haloes (Aq-B and Aq-E) have consistently steeper 
profiles and show no obvious break. Their densities in the Solar 
shell are none the less comparable to the many-progenitor haloes. 
Aq-F is dominated by a single progenitor, the debris of which re- 
tains a high degree of unmixed structure at z = (see also Fig. [15]). 

We show projected surface brightness profiles in Fig. [14] 
As with their three-dimensional counterparts, two characteristic 
shapes distinguish the many- and few-progenitor haloes. The few- 
progenitor haloes are centrally concentrated and well fit in their 
innermost ^ 10 kpc by Sersic profiles with 1.5 < n < 2.2. Be- 
yond 10 kpc, extended profiles with a more gradual rollover (de- 
scribed by Sersic profiles with n ^ 1 and 25 < res < 35 kpc) 



^ Both the Sagittarius and Milky Way halo stellar mass estimates are highly 
uncertain; it is unclear what contribution is made by the Sgr debris to esti- 
mates of the halo mass, although both the stream and the Virgo overdensity 
were masked out in the analysis of Bell et al. (2008) for which a value of 
~ 3 X 10^ Mq in the range 3 < r < 40 kpc was obtained from a broken 
power-law fit to the remaining 'smooth' halo. 
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Figure 13. Spherically averaged density profiles for our six stellar haloes 
(shown as thin lines below the = 7 radius of Navarro et al. 2008, at which 
the circular velocity of the dark matter halo has converged to an accuracy 
of 1%). Arrows mark the break radii of broken power-law fits to each pro- 
file. Dashed lines show Einasto profile fits to the corresponding dark matter 
haloes (Navarro et al. 2008). Grey vertical lines demarcate our outer halo 
region (dotted) and the Solar neighbourhood (solid); coloured vertical bars 
indicate r2oo for the dark haloes. For reference we overplot representative 
data for the Milky Way (orange): estimates of the halo density in the So- 
lar neighbourhood (symbols) from Gould, Flynn & Bahcall (1998, square) 
and Fuchs & JahreiB (1998, circle), and the best-fitting broken power-law 
of Bell et al. (excluding the Sagittarius stream and Virgo overdensity). 

are a better fit to the many-progenitor haloes. In their centres, how- 
ever, the many-progenitor haloes display a steep central inflection 
in surface brightness. As a consequence of these complex profiles, 
Sersic fits over the entire halo region (which we defined to begin 
at 3 kpc) are not fully representative in either case. To illustrate 
this broad dichotomy in Fig. O Sersic fits to a smoothly growing 
halo (Aq-C) beyond 10 kpc and a few-progenitor halo (Aq-E) in- 
terior to 10 kpc are shown. Abadi et al. (2006) found the average 
of their simulated stellar haloes to be well-fit by a Sersic profile 
(n = 6.3, reff = 7.7 kpc) in the radial range 30 < r < 130 kpc, 
which we show as an orange dashed line in Fig. [14] This profile is 
close to the 'mean' profile of our halos A, C and D interior to 30 kpc 
(neglecting the significant fluctuations and inflections within each 
individual halo in Fig.[T4]), but does not capture the sharp decline of 
our haloes at radii beyond 150 kpc. Fig.[T4]also shows (as dashed 
grey lines) the fits of Ibata et al. (2007) to the haloes of M31 (com- 
prising an r^/^ spheroid and shallow powerlaw tail at large radii) 
and M33 (powerlaw tail only). 

There is evidence for multiple kinematic and chemical subdi- 
visions within the Galactic globular cluster population (e.g. Searle 
& Zinn 1978; Frenk & White 1980; Zinn 1993; Mackey & Gilmore 
2004, and refs. therein). This has led to suggestions that at least 
some of these cluster subsets may have originated in accreted satel- 
lites (Bellazzini, Ferraro & Ibata 2003; Mackey & Gilmore 2004; 
Forbes, Strader, & Brodie 2004). Support for this conclusion in- 
cludes the presence of five globular clusters in the Fornax dwarf 
spheroidal (Hodge 1961) and the association of several Galac- 
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Figure 14. Radially averaged surface brightness profiles. Dashed lines show 
illustrative Sersic fits to haloes Aq-E and Aq-C (see text), with arrows in- 
dicating the corresponding scale radii. We show sections of equivalent pro- 
files for the haloes of M31 (including the inner r^/^ 'spheroid') and M33 
(beyond 10 kpc) as dashed grey lines (Ibata et al. 2007). We overplot the 
surface number density (right-hand axis) of globular clusters in M3 1 (yel- 
low squares) and the Milky Way (orange squares), with 40 and 10 clusters 
per bin, respectively. These profiles have been arbitrarily normalized to cor- 
respond to an estimate of the surface brightness of halo stars in the Solar 
neighbourhood from Morrison (1993), shown by a orange triangle. Vertical 
lines are as in Fig. [73] 

tic clusters with the Sagittarius nucleus and debris (e.g. Layden 
& Sarajedini 2000; Newberg et al. 2003; Bellazzini et al. 2003). 
Similarities with the 'structural' properties of stellar populations 
in the halo have motivated a longstanding interpretation of glob- 
ular clusters as halo (i.e. accretion debris) tracers (e.g. Lynden- 
Bell & Lynden-Bell 1995). We therefore plot in Fig. [14] the sur- 
face density profile of globular clusters in the Milky Way (Har- 
ris 1996) and M31 (confirmed GCs in the Revised Bologna Cat- 
alogue - RBC v3.5, March 2008 Galleti et al. 2004, 2006, 2007; 
Kim et al. 2007; Huxor et al. 2008). The Milky Way data have been 
projected along an arbitrary axis, and the normalization has been 
chosen to match the surface density of Milky Way clusters to an 
estimate of the surface brightness of halo stars in the Solar neigh- 
bourhood (/J.V = 27.7mag/arcsec^; Morrison 1993). We caution 
that the RBC incorporates data from ongoing surveys as it becomes 
available: the M31 GC profile shown here is therefore substantially 
incomplete, particularly with regard to the sky area covered beyond 
- 20-30 kpc. 

Abadi et al. (2006) showed that their average stellar halo Ser- 
sic fit also approximates the distribution of globular clusters in the 
Milky Way and M31. As stated above, the inner regions of our 
haloes Aq-A, Aq-C and Aq-D are in broad agreement with the 
Abadi et al. halo profile, and hence show some similarities with the 
observed globular cluster profiles also. Both the halo and cluster 
samples show strong variations from halo to halo, however, and the 
comparison of these small samples is inconclusive. A close corre- 
spondence between accreted halo stars and globular clusters would 
be expected only if the majority of clusters are accreted, if accreted 
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satellites contribute a number of clusters proportionate to their stel- 
lar mass, and if all stripped clusters have an equal probability of 
surviving to ^ = 0. None of these assumptions is realistic, and fur- 
ther work is required better to constrain the relationship between 
globular clusters and stellar haloes. 

The multicomponent nature of our haloes, which gives rise to 
the local structure in their overall profiles, is examined in more de- 
tail in Fig. [15] Here the density profiles of the major contributors 
shown in Fig.[TOlare plotted individually (progenitors contributing 
< 5% of the halo have been added to the panel for Aq-F). It is clear 
from these profiles that material from a given progenitor can be de- 
posited over a wide range of radii. The few-progenitor haloes show 
strong gradients in pr^ while more uniform distributions of this 
quantity are seen in their sub-dominant contributors and in most 
contributors to the many-progenitor haloes. 

Finally, we show in Fig.[T6lthe time at which the satellite pro- 
genitors of halo stars at a given radius were accreted (this infall 
time is distinct from the time at which the stars themselves were 
stripped, which may be considerably later). An analogous infall 
time can be defined for the surviving satellites, which are shown 
as points in Fig. [161 We would expect little information to be en- 
coded in an instantaneous sample of the radii of surviving satellites, 
but their infall times can none the less be usefully compared with 
those of halo stars. 

A gradient to earlier infall times with decreasing radius is ap- 
parent in both the satellites and the many-progenitor haloes. In the 
case of the haloes, this reflects the fact that relatively larger apocen- 
tres are associated with later-infalling satellites, which enable them 
to deposit material over a greater radial range. Assembly in this 
manner is arguably not adequately characterised as 'inside out' for- 
mation; late infalling material is added at all radii but has a greater 
maximum extent than earlier-infalling material. The result is that 
earlier-infalling material comes to dominate towards the centre. For 
the few-progenitor haloes the profile of infall time is essentially flat 
(or shows sharp transitions between populations), more closely re- 
flecting the contributions of individual progenitors. 

Further to our discussion of satellite survival in our haloes in 
Section |4l2l it is interesting that amongst the surviving satellites, we 
observe several accreted dX z > 1. For example, in the case of Aq- 
E, six surviving satellites are accreted dXz ^ 3.5; at the present day 
this group is found in association with a concentration of halo stars 
from a stellar halo progenitor also infalling at this time. The ma- 
jority of survivors in each halo are accreted recently, however, and 
typically more recently than the stellar halo progenitors. The op- 
posite is true for the earliest-accreted survivors, which are accreted 
earlier than the halo at the notably small radii at which they are now 
found. In general, at any given instant the majority of satellites are 
more likely to be located nearer to the apocentre of their orbit than 
the pericentre; furthermore, the orbits of the most massive satellites 
are likely to have been more circular than their disrupted siblings 
and dynamical friction may act to reinforce such a trend. There- 
fore, the locations of early-infalling survivors are likely to be fairly 
represented by their radius in Fig. [161 Dynamical friction acts to 
contract but also to circularize orbits. Plausibly these survivors are 
those that have sunk slowly as the result of their initially low orbital 
eccentricities. 



4.3.2 Stellar populations 

In this section, we show how the multicomponent nature of our 
stellar haloes is reflected in their metallicity profiles, and contrast 
the stellar populations of surviving satellites with those of halo 
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Figure 15. Individual density profiles (multiplied by r^) for stars con- 
tributed by each of the most significant progenitors of the halo (defined in 
Section (33). Line types indicate the rank order of a progenitor contribution: 
the bold coloured line in each panel indicates the most significant contribu- 
tor, while lesser contributions are shown by increasingly lighter and thinner 
lines. Vertical solid and dashed lines indicate the Solar shell and virial ra- 
dius respectively, as Fig.[T3l Individual stellar halo components contribute 
over a wide radial range, and different components 'dominate' at particular 
radii. This figure can be used to interpret the radial trends shown in other 
figures. 



progenitors. We caution that a full comparison of the relationship 
between the stellar halo and surviving satellites will require more 
sophisticated modelling of the chemical enrichment process than 
is included in our fiducial model, which adopts the instantaneous 
recycling approximation and does not follow individual elemen- 
tal abundances. We will address this detailed chemical modelling 
and related observational comparisons in a subsequent paper (De 
Lucia et al. in prep.). The model we adopt here tracks only total 
metallicity, defined as the total mass fraction of all metals relative 
to the Solar value, ^/Zq (the absolute value of which cannot be 
compared directly with measurements of [Fe/H]). This model can 
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Figure 16. Lines show, for halo stars at a given radius at z = 0, the mean 
(solid), median and 10/90th percentile (dotted) redshift at which their par- 
ent galaxy was accreted on to the main halo (not the time at which the 
stars themselves were stripped). Filled circles show the redshift at which 
surviving satellites were accreted; triangles indicate satellites accreted be- 
fore z = 7. Within the solar shell, the stellar halo is typically old in this 
'dynamical' sense, whereas beyond 100 kpc its young 'dynamical' age is 
comparable to that of the surviving satellite population. In many cases the 
innermost satellites represent a relic population that is 'older' than the stel- 
lar halo at comparable radii. 



nevertheless address the relative enrichment levels of different pop- 
ulations. 

Fig. [T7] shows the spherically averaged metallicity gradient 
in each halo. Our many-progenitor haloes are characterised by a 
metallicity distribution of width ^ 1 dex and approximately con- 
stant mean value, fluctuating by less than ±0.5 dex over a range of 
100 kpc. This is comparable to observations of the M31 halo, which 
show no significant gradient (metallicities varying by ±0.14 dex) 
in the range 30-60 kpc (Richardson et al. 2009). Localised struc- 
ture is most apparent in the few-progenitor haloes: Aq-F shows a 
clear separation into two components, while Aq-B and Aq-E ex- 



hibit global trends of outwardly declining metallicity gradients. In 
all cases the mean metallicity within the Solar radius is relatively 
high. These features can be explained by examining the relative 
weighting of contributions from individual progenitors at a given 
radius, as shown in the density profiles of Fig. [151 bearing in mind 
the mass-metallicity relation for satellites that arises in our model. 
Where massive progenitors make a significant luminosity- weighted 
contribution, the haloes are seen to be metal-rich. Overall, metallic- 
ity gradients are shallower in those haloes where many significant 
progenitors make a comparable contribution, smoothing the distri- 
bution over the extent of the halo. Conversely, metallicity gradients 
are steeper where only one or two disproportionately massive satel- 
lites make contributions to the halo (as indicated by the luminosity 
functions of Fig. [9]). Sharp contrasts are created between the radii 
over which this metal-rich material is deposited (massive satellites 
suffer stronger dynamical friction and sink more rapidly, favouring 
their concentration at the centres of haloes) and a background of 
metal-poor material from less massive halo progenitors. This ef- 
fect is clearly illustrated by the sharp transition in Aq-F and at two 
locations (centrally and at ^ 100 kpc) in Aq-E. 

It follows that the process by which our smooth haloes are as- 
sembled, which gives rise to the steep gradients of progenitor infall 
time with redshift shown in Fig. [T6j also acts to erase metallicity 
gradients. As a result, measurements of (for example) [Fe/H] alone 
do not constrain the local infall time; a metal-poor halo need not be 
'old' in the sense of early assembly. A particularly notable example 
of this is Aq-E, where the centrally dominant metal-rich material 
was assembled into the halo considerably earlier (z ^ 3) than the 
diffuse outer envelope of relatively metal-poor material (z ^ 1). 
This is a manifestation of a mass-metallicity relation in satellites: at 
fixed luminosity, an earlier infall time is 'compensated' for by more 
rapid star formation, resulting in a comparable degree of overall en- 
richment as that for a satellite with similar luminosity infalling at 
lower redshift. Abundance ratios such as [a/Fe] indicate the time 
taken by a given stellar population to reach its observed level of en- 
richment, and so distinguish between rapidly forming massive pop- 
ulations, truncated by early accretion to the halo, and populations 
reaching similar mass and metallicity through gradual star forma- 
tion (e.g. Shetrone, Cote & Sargent 2001; Tolstoy et al. 2003; Venn 
et al. 2004; Robertson et al. 2005). 

Fig. [18] shows luminosity-weighted metallicity distribution 
functions (MDFs) for two selections of halo stars: a 'Solar shell' 
(5 < r < 12 kpc; dashed lines) and the entire halo as defined in 
Section [331 (dotted). We compare these to MDFs for stars in the 
surviving satellites in each halo, separating bright (My < —10, 
r < 280 kpc; thick, coloured) and 'faint' (-10 < My < -5; 
thin, grey) subsets. All distributions are normalized individually to 
the total luminosity in their sample of stars. 

The MDF of Solar-shell halo stars is typically broad, and tends 
to peak at slightly higher metallicity (by < 0.5 dex) than the aggre- 
gated surviving bright satellites. The halo as a whole is compara- 
ble to the Solar shell. A clear disparity is only evident in Aq-E, 
where the halo appears to reflect more closely the distribution of 
fainter, lower-metallicity satellites. In all cases, the MDF of these 
faint satellites peaks at considerably lower metallicity than in the 
halo or brighter satellites. We find that the 'average' halo has an 
equivalent number of very metal-poor stars to the surviving bright 
satellites, although there are clear exceptions in individual cases. 
The fainter satellites have a substantially greater fraction of very 
metal-poor stars, in accordance with their low mean metallicities. 
Surviving satellites contain a greater fraction of moderately metal- 
poor stars (log^Q(Z/Z0) < —2.5) than the halo. 
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Figure 17. Radial profiles of luminosity- weighted metallicity (ratio of total Figure 18. Metallicity distribution functions of bright (My < -10; solid 

metal mass fraction to the Solar value) for spherical shells in our six haloes, coloured) and faint (-10 < My < -5; soHd grey) satellites, halo stars in 

showing the mean (solid) and median (thick dotted) profiles, bracketed by the 'Solar shell' (dashed) and the entire halo (3 < r < 280 kpc, dotted), 

the 10th and 90th percentiles (dotted). z is the total mass fraction of all metals. 



Our halo models suggest that similar numbers of comparably brightest satellites have a median age of 8 Gyr and a substantial 

luminous (and hence metal-rich) satellites contribute to the bright tail to young ages (with ^ 20% younger than 4 Gyr and ~ 90% 

end of both the halo-progenitor and the surviving- satellite luminos- younger than the median halo age). The distribution of old stars in 

ity functions, and that these bright satellites are the dominant con- the faintest surviving satellites is similar to that of the halo, 

tributors to the halo. This supports the view that halo MDFs should The true age distribution of halo stars is poorly constrained 

resemble those of bright survivor satellites in their metal-poor tails. in comparison to that of the satellites (e.g. Tolstoy, Hill & Tosi 

At very low metallicities the halo is dominated by the contribution 2009). By comparing the colour and metallicity distributions of 

of low-luminosity satellites which are exclusively metal-poor; the Milky Way halo stars to those of the Carina dSph, Unavane, Wyse 

stars associated with these faint contributors are expected to repre- & Gilmore (1996) have argued that similar satellites (i.e. those with 

sent only a very small fraction of the total halo luminosity. a substantial fraction of intermediate- age stars) could not contribute 

Finally, Fig. [19] compares the luminosity- weighted age distri- more than ^ 1% to the halo (equivalent to a maximum of ^ 60 halo 

butions of halo stars in the Solar shell with those in the surviving progenitors of Carina's luminosity). A corresponding limit of ^ 6 

satellites (My < —5), separated into bright and faint subsets. The Fornax-like accretions in the last ^ 10 Gyr was derived from an 

average of all six haloes contains essentially no stars younger than analysis of higher metallicity stars by the same authors, consistent 

5 Gyr (if we exclude halo Aq-F, which is strongly influenced by with the progenitor populations of our simulated stellar haloes, 

the late accretion of an SMC-like object, this minimum age rises to It is important in this context that the satellites themselves 

8 Gyr). The median age of halo stars is 11 Gyr. By contrast, the form hierarchically. In our models, between ten and twenty pro- 
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Figure 19. The cumulative luminosity-weighted age distribution (mean of 
all six simulations) for halo stars in the Solar shell (5 < r < 12 kpc, or- 
ange, top panel) compared to bright (—15 < My < —10; light green, 
bottom) and faint (—10 < My < —5; dark green, centre) satellites 
(My < —10), showing individual contributions from each halo (dashed, 
colours as in previous figures) to the mean value represented by each panel. 
The total stellar masses of these three components over all haloes are 
1.04 X 10^, 7.45 X 10^ and 3.45 x 10^ Mq, respectively 

genitors are typical for a (surviving) galaxy of stellar mass compa- 
rable to Sagittarius, or five to ten for a Fornax analogue. Satellites 
in this mass range are the most significant contributors to our stellar 
haloes. Their composite nature is likely to be reflected in their stel- 
lar population mix and physical structure, which could complicate 
attempts to understand the halo 'building blocks' and the surviving 
satellites in terms of simple relationships between mass, age and 
metallicity. 



5 CONCLUSIONS 

We have presented a technique for extracting information on the 
spatial and kinematic properties of galactic stellar haloes that com- 
bines a very high resolution fully cosmological ACDM simulation 
with a semi-analytic model of galaxy formation. We have applied 
this technique to six simulations of isolated dark matter haloes sim- 
ilar to or slightly less massive than that of the Milky Way, adopt- 
ing a fiducial set of paramter values in the semi-analytic model 
GALFORM. The structural properties of the surviving satellites have 
been used as a constraint on the assignment of stellar populations 
to dark matter. We found that this technique results in satellite pop- 
ulations and stellar haloes in broad agreement with observations of 
the Milky Way and M31, if allowance is made for differences in 
dark halo mass. 

Our method of assigning stellar populations to dark matter 
particles is, of course, a highly simplified approach to modelling 
star formation and stellar dynamics. The nature of star formation 
in dwarf galaxy haloes remains largely uncertain. In future, ob- 
servations of satellites interpreted alongside high-resolution hydro- 
dynamical simulations will test the validity of approaches such as 



ours. As a further simplification, our models do not account for 
a likely additional contribution to the halo from scattered in situ 
(disc) stars, although we expect this contribution to be minimal far 
from the bulge and the disc plane. The results outlined here there- 
fore address the history, structure and stellar populations of the ac- 
creted halo component in isolation. 

Our results can be summarised as follows: 

• Our six stellar haloes are predominantly built by satellite ac- 
cretion events occurring between 1 < 2; < 3. They span a range 
of assembly histories, from 'smooth' growth (with a number of 
roughly equally massive progenitors accreted steadily over a Hub- 
ble time) to growth in one or two discrete events. 

• Stellar haloes in our model are typically built from fewer 
than 5 significant contributors. These significant objects have stel- 
lar masses comparable to the brightest classical dwarf spheroidals 
of the Milky Way; by contrast, fewer faint satellites contribute to 
the halo than are present in the surviving population. 

• Typically, the most massive halo contributor is accreted at a 
lookback time of between 7 and 1 1 Gyr (z ^ 1.5 — 3) and deposits 
tidal debris over a wide radial range, dominating the contribution at 
large radii. Stars stripped from progenitors accreted at even earlier 
times usually dominate closer to the centre of the halo. 

• A significant fraction of the stellar halo consists of stars 
stripped from individual surviving galaxies, contrary to expecta- 
tions from previous studies (e.g. Bullock & Johnston 2005). It is 
the most recent (and significant) contributors that are likely to be 
identifiable as surviving bound cores. Such objects have typically 
lost ^ 90% of their original stellar mass. 

• We find approximately power-law density profiles for the stel- 
lar haloes in the range 10 < r < 100 kpc. Those haloes formed 
by a superposition of several comparably massive progenitors have 
slopes similar to those suggested for the Milky Way and M31 
haloes, while those dominated by a disproportionally massive pro- 
genitor have steeper slopes. 

• Our haloes have strongly prolate distributions of stellar mass 
in their inner regions (c/a ^ 0.3), with one exception, where an 
oblate, disc-like structure dominates the inner 10 — 20 kpc. 

• Haloes with several significant progenitors show little or no 
radial variation in their mean metallicity (Z/Zq) up to 200 kpc. 
Those in which a small number of progenitors dominate show 
stronger metallicity gradients over their full extent or sharp transi- 
tions between regions of different metallicity. The centres of these 
haloes are typically more enriched than their outer regions. 

• The stellar populations of the halo are likely to be chemically 
enriched to a level comparable to that of the bright surviving satel- 
lites, but to be as old as the more metal-poor surviving 'ultra faint' 
galaxies. The very metal-poor tail of the halo distribution is dom- 
inated by contributions from a plethora of faint galaxies that are 
insignificant contributors to the halo overall. 
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